drop _all
clear
set more off
global DFAD		"/Mac/DFAD"
global input 	"$DFAD/DFAD-INPUT"
global output 	"$DFAD/DFAD-OUTPUT"
global paper 	"$DFAD/DFAD-PAPER/Latex"
cd "$DFAD"

*** Purely cosmetic -- sets global macros that define colors in graphs below
	do "/Mac/Macros/JCS Graph Colors.do"
	mycolors white

************
***	DTP3 ***
************

	use "$output/JDS/dtp3.dta", clear
	
*** Changes
	foreach var in year dhs who {
		sort country year
		by country: gen lag`var' = `var'[_n-1]
		by country: gen fwd`var' = `var'[_n+1]
	}
	foreach var in dhs who{
		gen D`var' = (`var' - lag`var')
		gen D10`var' = 10*(`var' - lag`var')/(year - lagyear)
		gen gr_`var' = 100*((`var'/lag`var')^(1/(year-lagyear)) - 1)
		gen fgr_`var' = 100*((fwd`var'/`var')^(1/(fwdyear-year)) - 1)
	}
	
*** Compare
	tostring year, g(mylabel)
	tostring lagyear, g(slagyear)
	replace wb_ccode = wb_ccode + ", " 
	gen		newcode = wb_ccode + mylabel
	replace mylabel = slagyear + "-" + mylabel
	gen 	type = 1
	replace type = 2 if Dwho>Ddhs 
	gen disease = "DTP3"
	gen ratio = who/dhs
	tempfile dtp
	save `dtp', replace
	
	replace country = "DRC" if country=="Congo Democratic Republic" | country=="Democratic Republic of the Congo"
	
	keep if africa==1
	reg ratio year if year<2000
	set obs 182
	replace year = 1999.9 if _n==182
	predict rhat
	replace ratio = rhat if _n==182

	twoway 	(lfitci ratio year if year<2000,   lcolor($c1) acolor(gs12) fcolor(gs12))									///
			(lfitci ratio year if year>=2000,  lcolor($c1) acolor(gs12) fcolor(gs12))									///
			(scatter ratio year if year<2000,  mcolor($c1) msymbol(Oh))													///
			(scatter ratio year if year>=2000, mcolor($c1) msymbol(X))													///
			(scatter ratio year if year<2000 & (ratio>1.2 | ratio<.9),  mcolor($c1) mlabcolor($c1) mlabel(country) msymbol(Oh))		///
			(scatter ratio year if year>=2000 & (ratio>1.2 | ratio<.9), mcolor($c1) mlabcolor($c1) mlabel(country) mlabpos(9) msymbol(X)),		///
			xlabel(, $myxlabel) 																						///
			ylabel(0.6 0.8 1.0 1.2 1.4 1.6 1.8, $myylabel) 																///
			xscale(range(2011) lcolor($cf))																				///
			yscale(range(0.6 1.8) lcolor($cf))																			///
			xtitle("", color($cf)) 																						///
			$myregion 																									///
			ytitle("Ratio of WHO to DHS coverage", color($cf))															///
			legend(off) 				
			*graph export "$paper/tsline_DTP_new.pdf", as(pdf) replace   
	
	
***************
***	MEASLES ***
***************

	use "$output/JDS/measles.dta", clear

*** Changes
	foreach var in year dhs who {
		sort country year
		by country: gen lag`var' = `var'[_n-1]
		by country: gen fwd`var' = `var'[_n+1]
	}
	foreach var in dhs who {
		gen D`var' = (`var' - lag`var')
		gen D10`var' = 10*(`var' - lag`var')/(year - lagyear)
		gen gr_`var' = 100*((`var'/lag`var')^(1/(year-lagyear)) - 1)
		gen fgr_`var' = 100*((fwd`var'/`var')^(1/(fwdyear-year)) - 1)
	}

*** Compare
	tostring year, g(mylabel)
	tostring lagyear, g(slagyear)
	replace wb_ccode = wb_ccode + ", " 
	gen	newcode = wb_ccode + mylabel
	replace mylabel = slagyear + "-" + mylabel
	gen 	type = 1
	replace type = 2 if Dwho>Ddhs 
	gen disease = "Measles"
	gen ratio = who/dhs
	tempfile measles
	save `measles', replace

	replace country = "DRC" if country=="Congo Democratic Republic" | country=="Democratic Republic of the Congo"

	keep if africa==1
	reg ratio year if year<2000
	set obs 182
	replace year = 1999.9 if _n==182
	predict rhat
	replace ratio = rhat if _n==182

	twoway 	(lfitci ratio year if year<2000,   lcolor($c1) acolor(gs12) fcolor(gs12))									///
			(lfitci ratio year if year>=2000,  lcolor($c1) acolor(gs12) fcolor(gs12))									///
			(scatter ratio year if year<2000,  mcolor($c1) msymbol(Oh))													///
			(scatter ratio year if year>=2000, mcolor($c1) msymbol(X))													///
			(scatter ratio year if year<2000 & (ratio>1.2|ratio<.8),  mcolor($c1) mlabcolor($c1) mlabel(country) msymbol(Oh))		///
			(scatter ratio year if year>=2000 & (ratio>1.2|ratio<.8), mcolor($c1) mlabcolor($c1) mlabel(country) mlabpos(9) msymbol(X)),		///
			xlabel(, $myxlabel) 																						///
			ylabel(0.6 0.8 1.0 1.2 1.4 1.6 1.8, $myylabel) 																///
			xscale(range(2011) lcolor($cf))																				///
			yscale(range(0.6 1.8) lcolor($cf))																			///
			xtitle("", color($cf)) 																						///
			$myregion 																									///
			ytitle("Ratio of WHO to DHS coverage", color($cf))															///
			legend(off) 				
			*graph export "$paper/tsline_measles_new.pdf", as(pdf) replace   
	 

********************************************
*** Regressions using both DTP & measles ***
********************************************

		use 	`dtp', clear
		append 	using 	`measles'
		gen 	dtp 		= disease=="DTP3"
		gen 	dhsdtp 		= dhs*dtp
		gen 	y2k 		= year>=2000
		gen 	dtpy2k		= dtp*y2k
		gen 	Lim 		= year>2008
		gen 	dtpLim 		= dtp*Lim
		label 	variable 	dhs 	"DHS rate"
		label 	variable 	dtp 	"DTP3"
		label 	variable 	dhsdtp 	"DHS rate $ \times $ DTP3"
		label 	variable 	y2k 	"Year $ \geq $ 2000"
		label 	variable 	dtpy2k 	"(Year $ \geq $ 2000) $ \times $ DTP3"
		label 	variable 	Lim		"Year $ \geq $ 2008"
		label 	variable 	dtpLim 	"(Year $ \geq $ 2008) $ \times $ DTP3"
	
		gen 	dhslevel	= dhs
		encode 	country, g(ccode)
		
		keep if africa==1
	

	*** Regressions
		reg who  y2k Lim		if dtp==1, cluster(ccode)
			est2vec health_140503, replace vars(y2k Lim dhs dtp dtpy2k dtpLim _cons) name(DTP1)
		reg who  y2k Lim dhs 	if dtp==1, cluster(ccode)
			est2vec health_140503, addto(health_140503) name(DTP2)
			replace dhs = Ddhs
		reg Dwho y2k Lim dhs 	if dtp==1, cluster(ccode)
			est2vec health_140503, addto(health_140503) name(DTP3)
			replace dhs = dhslevel
		reg who  y2k Lim		if dtp==0, cluster(ccode)
			est2vec health_140503, addto(health_140503) name(M1)
		reg who  y2k Lim dhs 	if dtp==0, cluster(ccode)
			est2vec health_140503, addto(health_140503) name(M2)
			replace dhs = Ddhs
		reg Dwho y2k Lim dhs 	if dtp==0, cluster(ccode)
			est2vec health_140503, addto(health_140503) name(M3)
			replace dhs = dhslevel
		reg who  y2k Lim 		dtp	dtpy2k dtpLim, cluster(ccode)	
			est2vec health_140503, addto(health_140503) name(BOTH1)
		reg who  y2k Lim dhs 	dtp dtpy2k dtpLim, cluster(ccode)
			est2vec health_140503, addto(health_140503) name(BOTH2)
			replace dhs = Ddhs
		reg Dwho y2k Lim dhs 	dtp dtpy2k dtpLim, cluster(ccode)
			est2vec health_140503, addto(health_140503) name(BOTH3)
			est2tex health_140503, preserve path("$paper") levels (90 95 99) extracol(2 3 5 6 8) fancy label replace digits(2) mark(stars)

exit	
	
